import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
Basemap is a matplotlib extension used to visualize and create geographical maps in python.
If you need further information on basemap, please visit basemap page: https://matplotlib.org/basemap/index.html
conda install -c anaconda basemap
# Import the Basemap and Matplotlib libraries:
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
Coastlines
First, let’s initialize a map with Basemap() function
Then, use the drawcoastlines() function to add coastlines on the map
fig = plt.figure(figsize = (8,8))
m = Basemap()
m.drawcoastlines()
plt.title("Coastlines", fontsize=20)
plt.show()
drawcoastlines function has the following main arguments:
linewidth: 1.0, 2.0, 3.0…
linestyle: solid, dashed…
color: black, red…
fig = plt.figure(figsize = (6,6))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='dashed', color='red')
plt.title("Coastlines", fontsize=20)
plt.show()
Countries
drawcountries() function to add countries on the mapfig = plt.figure(figsize = (8,8))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries()
plt.title("Country boundaries", fontsize=20)
plt.show()
fig = plt.figure(figsize = (12,12))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries(linewidth=1.0, linestyle='solid', color='r')
plt.title("Country boundaries", fontsize=20)
plt.show()
Draw major rivers
Use the drawrivers() function to add major rivers on the map.
drawrivers() function can take linewidth, linestyle, and color arguments
fig = plt.figure(figsize = (12,12))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries(linewidth=1.0, linestyle='solid', color='k')
m.drawrivers(linewidth=0.5, linestyle='solid', color='#0000ff')
plt.title("Major rivers", fontsize=20)
plt.show()
Fill continents
This function is used to draw color filled continents
Use fillcontinents() function to fill continents
fig = plt.figure(figsize = (8,8))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries(linewidth=1.0, linestyle='solid', color='k')
m.fillcontinents()
plt.title("Color filled continents", fontsize=20)
plt.show()
The fillcontinents() function can take the following arguments:
color: fills continents (default gray)
lake_color: fills inland lakes
alpha: sets transparency for continents
fig = plt.figure(figsize = (12,12))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries(linewidth=1.0, linestyle='solid', color='k')
m.fillcontinents(color='coral',lake_color='aqua', alpha=0.9)
plt.title("Color filled continents", fontsize=20)
plt.show()
Draw map boundary
The drawmapboundary() function is used to draw the earth boundary on the map.
drawmapboundary() function can take the following arguments:
linewidth: sets line width for boundary line (default: 1)
color: sets the color of the boundary line (default: black)
fill_color: fills the map background region
fig = plt.figure(figsize = (6,6))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries(linewidth=1.0, linestyle='solid', color='k')
m.fillcontinents(color='coral',lake_color='aqua')
m.drawmapboundary(color='b', linewidth=2.0, fill_color='aqua')
plt.title("Filled map boundary", fontsize=20)
plt.show()
Draw and label longitude lines
drawmeridians() function is used to draw & label meridians/longitude lines
drawmeridians()function can take the following arguments:
range() for integer values & np.arange() for floats values.color: sets the color of the line longitude lines (default: black).
textcolor: sets the color of labels (default: black)
linewidth: sets the line width for the longitude lines
dashes: sets the dash pattern for the longitude lines (default: [1,1])
labels: sets the label’s position with four values [0,0,0,0] representing left, right, top, & bottom. Change these values to 1 where you want the labels to appear.
fig = plt.figure(figsize = (8,8))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries(linewidth=1.0, linestyle='solid', color='k')
m.fillcontinents(color='coral',lake_color='aqua')
m.drawmeridians(range(0, 360, 30), color='k', linewidth=1.0, dashes=[4, 4], labels=[1, 0, 0, 1])
plt.title("Longitude lines", fontsize=20, pad=30)
plt.show()
Draw and label latitude lines
drawparallels() function is used to draw & label parallels/latitude lines.
drawparallels()function uses similar arguments like drawmeridians()
fig = plt.figure(figsize = (8,8))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries(linewidth=1.0, linestyle='solid', color='k')
m.fillcontinents(color='coral',lake_color='aqua')
m.drawparallels(range(-90, 100, 20), color='k', linewidth=1.0, dashes=[4, 4], labels=[1, 1, 0, 0])
plt.title("Latitude lines", fontsize=20)
plt.show()
Let's put the drawmeridians and drawparallels functions together
fig = plt.figure(figsize = (12,12))
m = Basemap()
m.drawcoastlines(linewidth=1.0, linestyle='solid', color='black')
m.drawcountries(linewidth=1.0, linestyle='solid', color='k')
m.fillcontinents(color='coral',lake_color='aqua')
m.drawmeridians(range(0, 360, 20), color='k', linewidth=1.0, dashes=[4, 4], labels=[0, 0, 0, 1])
m.drawparallels(range(-90, 100, 10), color='k', linewidth=1.0, dashes=[4, 4], labels=[1, 0, 0, 0])
plt.ylabel("Latitude", fontsize=15, labelpad=35)
plt.xlabel("Longitude", fontsize=15, labelpad=20)
plt.show()
Projection, bounding box, & resolution
projection, bounding box, & resolution of a mapMap projection:
Inside the Basemap() function, the projection=" argument can take several pre-defined projections
To specify the desired projection, use the general syntax shown below:
m = Basemap(projection='aeqd')
m = Basemap(projection='cyl')
| Projection name | Basemap syntax | | Projection name | Basemap syntax | |----------------------------------- |---------------- |--- |----------------------------------- |---------------- | | Equidistant Conic | eqdc | | Cylindrical Equal Area | cea | | Miller Cylindrical | mill | | McBryde-Thomas Flat-Polar Quartic | mbtfpq | | Mercator | merc | | Azimuthal Equidistant | aeqd | | Stereographic | stere | | Sinusoidal | sinu | | Polyconic | poly | | Rotated Pole | rotpole | | Oblique Mercator | omerc | | Cylindrical Equidistant | cyl | | Gnomonic | gnom | | North-Polar Stereographic | npstere | | Mollweide | moll | | South-Polar Stereographic | spstere | | Lambert Conformal | lcc | | Hammer | hammer | | Transverse Mercator | tmerc | | Geostationary | geos | | North-Polar Lambert Azimuthal | nplaea | | Near-Sided Perspective | nsper | | Gall Stereographic Cylindrical | gall | | Eckert IV | eck4 | | North-Polar Azimuthal Equidistant | npaeqd | | Albers Equal Area | aea | | Kavrayskiy VII | kav7 | | South-Polar Azimuthal Equidistant | spaeqd | | Orthographic | ortho | | Cassini-Soldner | cass | | van der Grinten | vandg | | Lambert Azimuthal Equal Area | laea | | South-Polar Lambert Azimuthal | splaea | | Robinson | robin |
Some projections require setting bounding box, map center, & map size of the map using the following arguments:
a) Bounding box/map corners:
llcrnrlat: lower-left corner geographical latitude
urcrnrlat: upper-right corner geographical latitude
llcrnrlon: lower-left corner geographical longitude
urcrnrlon: upper-right corner geographical longitude
Example:
m = Basemap(projection='cyl', ,llcrnrlat=-80,urcrnrlat=80,llcrnrlon=-180,urcrnrlon=180)
b) Map center:
lon_0: central longitude
lat_0: central latitude
Example:
m = Basemap(projection='ortho', lon_0 = 25, lat_0 = 10)
c) Map resolution: The map resolution argument determines the quality of vector layers such as coastlines, lakes, & rivers etc.
The available options are:
pip install basemap-data-hires
conda install -c conda-forge basemap-data-hires
# Create a global map with a Mercator Projection
fig = plt.figure(figsize = (8,8))
m = Basemap(projection='merc',llcrnrlat=-80,urcrnrlat=80,llcrnrlon=-180,urcrnrlon=180)
m.drawcoastlines()
m.fillcontinents(color='tan',lake_color='lightblue')
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawmapboundary(fill_color='lightblue')
plt.title("Mercator Projection", fontsize=20)
Text(0.5, 1.0, 'Mercator Projection')
# Create a global map with a Cylindrical Equidistant Projection
fig = plt.figure(figsize = (8,8))
m = Basemap(projection='cyl',llcrnrlat=-80,urcrnrlat=80,llcrnrlon=-180,urcrnrlon=180)
m.drawcoastlines()
m.fillcontinents(color='tan',lake_color='lightblue')
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawmapboundary(fill_color='lightblue')
plt.title(" Cylindrical Equidistant Projection", fontsize=20)
Text(0.5, 1.0, ' Cylindrical Equidistant Projection')
# Create a global map with Orthographic Projection
fig = plt.figure(figsize = (4,4))
m = Basemap(projection='ortho', lon_0 = 25, lat_0 = 10)
m.drawcoastlines()
m.fillcontinents(color='tan',lake_color='lightblue')
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawmapboundary(fill_color='lightblue')
plt.title("Orthographic Projection", fontsize=18)
Text(0.5, 1.0, 'Orthographic Projection')
# Create a global map with a Robinson Projection
fig = plt.figure(figsize = (8,8))
m = Basemap(projection='robin',llcrnrlat=-80,urcrnrlat=80,llcrnrlon=-180,urcrnrlon=180, lon_0 = 0, lat_0 = 0)
m.drawcoastlines()
m.fillcontinents(color='tan',lake_color='lightblue')
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawmapboundary(fill_color='lightblue')
plt.title(" Robinson Projection", fontsize=20)
Text(0.5, 1.0, ' Robinson Projection')
fig = plt.figure(figsize = (8,8))
m = Basemap(projection='cyl',llcrnrlon=32.5,llcrnrlat=3,urcrnrlon=49,urcrnrlat=15, resolution = 'l')
m.drawcoastlines()
m.fillcontinents(color='tan',lake_color='lightblue')
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawmapboundary(fill_color='lightblue')
m.drawcoastlines()
plt.show()
b) By passing the central longitude & central latitude, as well as the width, & height values.
fig = plt.figure(num=None, figsize=(8, 8) )
m = Basemap(projection='aea', width=1700000, height=1500000, resolution='l',lon_0= 40.5,lat_0=8.7)
m.fillcontinents(color='tan',lake_color='aqua')
m.drawmapboundary(fill_color='lightblue')
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawcoastlines()
plt.show()
Basemap provides several background relief images including:
a) Land-sea mask image
b) NASA blue marble
c) NOAA etopo
d) Shaded relief
drawlsmask() function used to draw lakes, land, & oceans at once.
It avoids fillcontinents & drawmapboundary functions.
drawlsmask() function can take the following arguments:
land_color: sets the color of the land (by default it’s gray)
ocean_color: sets the color of the oceans (by default it’s white)
resolution: sets the coastline resolution (by default it’s l, you can change it to: c, l, i , h, or f )
lakes: plots lakes & ponds (by it default it’s True)
grid: sets the land/sea mask grid spacing in minutes (by default it’s 5; You can change it to: 10, 2.5, & 1.25)
fig = plt.figure(figsize = (8,8))
m = Basemap(projection='cyl',llcrnrlon=32.5,llcrnrlat=3,urcrnrlon=49,urcrnrlat=15, resolution = 'i')
m.drawlsmask(land_color = "#ddaa66", ocean_color="#7777ff", resolution = 'i', lakes=True, grid=1.25)
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
plt.title("Land-sea mask image", fontsize=20)
plt.show()
Use the bluemarble() argument to apply the NASA blue marble image as a background map
To downgrade or upgrade image resolution use the scale argument
fig = plt.figure(figsize = (8,8))
m = Basemap(projection='cyl',llcrnrlon=32.5,llcrnrlat=3,urcrnrlon=49,urcrnrlat=15, resolution='i')
m.bluemarble(scale=1.0)
m.drawcoastlines()
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
plt.title("NASA Blue Marble image as background map", fontsize=18)
Text(0.5, 1.0, 'NASA Blue Marble image as background map')
Use shadedrelief() argument to apply shaded relief image as a background map
fig = plt.figure(figsize = (8,8))
m = Basemap(projection='cyl',llcrnrlon=32.5,llcrnrlat=3,urcrnrlon=49,urcrnrlat=15, resolution='i')
m.shadedrelief()
m.drawcoastlines()
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
plt.title("Shaded relief image as background map", fontsize=18)
Text(0.5, 1.0, 'Shaded relief image as background map')
Use etopo() argument to apply shaded relief image as a background map
fig = plt.figure(figsize = (8,8))
m = Basemap(projection='cyl',llcrnrlon=32.5,llcrnrlat=3,urcrnrlon=49,urcrnrlat=15, resolution='i')
m.etopo(scale=1.2)
m.drawcoastlines()
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
plt.title("Shaded relief image as background map", fontsize=18)
Text(0.5, 1.0, 'Shaded relief image as background map')
To read the netcdf data, we’ll use the Dataset class from the netcdf4-python library
from netCDF4 import Dataset as dataset
nc = dataset('ECMWF_temp2m.nc')
nc.dimensions
{'longitude': <class 'netCDF4._netCDF4.Dimension'>: name = 'longitude', size = 144,
'latitude': <class 'netCDF4._netCDF4.Dimension'>: name = 'latitude', size = 73,
'time': <class 'netCDF4._netCDF4.Dimension'>: name = 'time', size = 62}
nc.variables
{'longitude': <class 'netCDF4._netCDF4.Variable'>
float32 longitude(longitude)
_FillValue: nan
units: degrees_east
long_name: longitude
unlimited dimensions:
current shape = (144,)
filling on,
'latitude': <class 'netCDF4._netCDF4.Variable'>
float32 latitude(latitude)
_FillValue: nan
units: degrees_north
long_name: latitude
unlimited dimensions:
current shape = (73,)
filling on,
'time': <class 'netCDF4._netCDF4.Variable'>
int32 time(time)
long_name: time
units: hours since 1900-01-01
calendar: proleptic_gregorian
unlimited dimensions:
current shape = (62,)
filling on, default _FillValue of -2147483647 used,
'p2t': <class 'netCDF4._netCDF4.Variable'>
int16 p2t(time, latitude, longitude)
_FillValue: -32767
units: K
long_name: 2 metre temperature
add_offset: 262.39847874753474
scale_factor: 0.0018355835199370632
missing_value: -32767
unlimited dimensions:
current shape = (62, 73, 144)
filling on}
lat = nc.variables['latitude'][:]
lon = nc.variables['longitude'][:]
time = nc.variables['time'][:]
t2 = nc.variables['p2t'][:]
Use the contourf() function to draw filled contour maps
Create a 2D meshgridmatrix from 1D longitude & latitude arrays
contourf() function mainly takes:
fig = plt.figure(num=None, figsize=(8, 8))
m = Basemap(projection='cyl', llcrnrlon=32.5, llcrnrlat=3, urcrnrlon=49, urcrnrlat=15, resolution='i')
x, y = m(*np.meshgrid(lon,lat))
cs = m.contourf(x, y ,np.squeeze(t2[0,:,:]), levels = 50, cmap=plt.cm.jet)
m.drawcoastlines()
m.drawmapboundary()
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawmeridians(range(33, 48, 2), color='k', linewidth=1.0, dashes=[4, 4], labels=[0, 0, 0, 1])
m.drawparallels(range(3, 15, 2), color='k', linewidth=1.0, dashes=[4, 4], labels=[1, 0, 0, 0])
plt.ylabel("Latitude", fontsize=15, labelpad=35)
plt.xlabel("Longitude", fontsize=15, labelpad=20)
cbar = m.colorbar(cs, location='right', pad="3%")
cbar.set_label('Temperature (K)', fontsize=13)
plt.title('2-meter temperature filled contour map', fontsize=15)
plt.show()
Use the contour() function to draw contour maps.
contour() function uses similar arguments like contourf()
fig = plt.figure(num=None, figsize=(8, 8))
m = Basemap(projection='cyl', llcrnrlon=32.5, llcrnrlat=3, urcrnrlon=49, urcrnrlat=15, resolution='i')
x, y = m(*np.meshgrid(lon,lat))
cs = m.contour(x, y ,np.squeeze(t2[4,:,:]), levels = 25 ,cmap=plt.cm.jet)
m.drawcoastlines()
m.drawmapboundary()
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawmeridians(range(33, 48, 2), color='k', linewidth=1.0, dashes=[4, 4], labels=[0, 0, 0, 1])
m.drawparallels(range(3, 15, 2), color='k', linewidth=1.0, dashes=[4, 4], labels=[1, 0, 0, 0])
plt.ylabel("Latitude", fontsize=15, labelpad=35)
plt.xlabel("Longitude", fontsize=15, labelpad=20)
cbar = m.colorbar(cs, location='right', pad="3%")
cbar.set_label('Temperature (K)', fontsize=13)
plt.clabel(cs, inline=True, fmt='%1.0f', fontsize=12, colors='k')
plt.title('2-meter temperature contour map', fontsize=15)
plt.show()
fig = plt.figure(num=None, figsize=(8, 8))
m = Basemap(projection='cyl', llcrnrlon=32.5, llcrnrlat=3, urcrnrlon=49, urcrnrlat=15, resolution='i')
x, y = m(*np.meshgrid(lon,lat))
clevs = np.arange(295, 321, 3)
#clevs = [295, 297, 299, 301, 303, 305, 307, 309, 311, 3315, 320]
cs = m.contour(x, y ,np.squeeze(t2[4,:,:]), clevs ,cmap=plt.cm.jet)
m.drawcoastlines()
m.drawmapboundary()
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
m.drawmeridians(range(33, 48, 2), color='k', linewidth=1.0, dashes=[4, 4], labels=[0, 0, 0, 1])
m.drawparallels(range(3, 15, 2), color='k', linewidth=1.0, dashes=[4, 4], labels=[1, 0, 0, 0])
plt.ylabel("Latitude", fontsize=15, labelpad=35)
plt.xlabel("Longitude", fontsize=15, labelpad=20)
cbar = m.colorbar(cs, location='right', pad="3%")
cbar.set_label('Temperature (K)', fontsize=13)
plt.clabel(cs, inline=True, fmt='%1.0f', fontsize=12, colors='k')
plt.title('2-meter temperature contour map', fontsize=15)
plt.show()
Use the pcolor() or pcolormesh() function to draw pseudo color maps
fig = plt.figure(figsize=(8, 8) )
m = Basemap(projection='cyl', llcrnrlon=32.5, llcrnrlat=3, urcrnrlon=49, urcrnrlat=15, resolution='i')
# cs = m.pcolor(x, y ,np.squeeze(t2[4,:,:]))
cs = m.pcolormesh(x, y ,np.squeeze(t2[4,:,:]))
# pcolor(): draw a pseudocolor plot.
# pcolormesh(): draw a pseudocolor plot (faster version for regular meshes).
m.drawmeridians(range(33, 48, 2), color='k', linewidth=1.0, dashes=[4, 4], labels=[0, 0, 0, 1])
m.drawparallels(range(3, 15, 2), color='k', linewidth=1.0, dashes=[4, 4], labels=[1, 0, 0, 0])
m.drawcoastlines()
m.drawmapboundary()
m.drawcountries(linewidth=1, linestyle='solid', color='k' )
plt.ylabel("Latitude", fontsize=15, labelpad=35)
plt.xlabel("Longitude", fontsize=15, labelpad=20)
cbar = m.colorbar(cs, location='right', pad="3%")
cbar.set_label('Temperature (K)', fontsize=13)
plt.title('2-meter temperature pseudo-color map', fontsize=15)
plt.show()
# Define the data points
latitudes = [10.90, 3.82, 6.60, 12.68, 8.09, 12.51, 14.21, 11.38, 4.19, 13.33]
longitudes = [34.41, 36.36, 40.55, 47.89, 44.39, 44.13, 44.9, 40.36, 39.14, 47.05]
data = [200, 150, 175, 189, 100, 315, 222, 275, 155, 190]
fig = plt.figure(figsize=(8, 8) )
# Create a Basemap instance
m = Basemap(
llcrnrlon=33, llcrnrlat=3, urcrnrlon=48, urcrnrlat=15,
projection='merc', resolution='i')
# Draw coastlines, countries, and states
m.drawcoastlines(linewidth=0.5)
m.drawcountries(linewidth=0.5)
m.drawstates(linewidth=0.5)
# Draw parallels and meridians
m.drawparallels(np.arange(0, 90, 4), labels=[1,0,0,0], fontsize=10)
m.drawmeridians(np.arange(-180, 180, 4), labels=[0,0,0,1], fontsize=10)
# Convert latitudes and longitudes to map coordinates
x, y = m(longitudes, latitudes)
# Plot the data points as markers on the map
m.scatter(x, y, c=data, cmap='cool', s=80, edgecolors='k', linewidths=0.5)
# Add a colorbar
plt.colorbar(label='Data', shrink =0.5)
# Add a title
plt.title('Data Points')
# Show the plot
plt.show()
fig = plt.figure(figsize=(8, 8) )
# Define the data points
latitudes = [10.90, 3.82, 6.60, 12.68, 8.09, 12.51, 14.21, 11.38, 4.19, 13.33]
longitudes = [34.41, 36.36, 40.55, 47.89, 44.39, 44.13, 44.9, 40.36, 39.14, 47.05]
data = [200, 150, 175, 189, 100, 315, 222, 275, 155, 190]
# Create a Basemap instance
m = Basemap(
llcrnrlon=33, llcrnrlat=3, urcrnrlon=48, urcrnrlat=15,
projection='merc', resolution='i'
)
# Draw coastlines, countries, and states
m.drawcoastlines(linewidth=0.5)
m.drawcountries(linewidth=0.5)
m.drawstates(linewidth=0.5)
# Draw parallels and meridians
m.drawparallels(np.arange(0, 90, 5), labels=[1,0,0,0], fontsize=10)
m.drawmeridians(np.arange(-180, 180, 5), labels=[0,0,0,1], fontsize=10)
# Convert latitudes and longitudes to map coordinates
x, y = m(longitudes, latitudes)
# Plot the data points as markers on the map
m.scatter(x, y, c=data, cmap='cool', s=100, edgecolors='k', linewidths=0.5)
# Add labels to the data points
for i in range(len(latitudes)):
plt.text(x[i]+50000, y[i]+50000, str(data[i]), fontsize=10, fontweight='bold', ha='center', va='center')
# Add a colorbar
plt.colorbar(label='Data', shrink = 0.5)
# Add a title
plt.title('Data Points Labled ')
# Show the plot
plt.show()